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Piecewise smooth hybrid systems, involving continuous and discrete variables, are suitable models 
for describing the multiscale regulatory machinery of the biological cells. In hybrid models, the 
discrete variables can switch on and off some molecular interactions, simulating cell progression 
through a series of functioning modes. The advancement through the cell cycle is the archetype 
of such an organized sequence of events. We present an approach, inspired from tropical geome- 
try ideas, allowing to reduce, hybridize and analyse cell cycle models consisting of polynomial or 
rational ordinary differential equations. 

1 Introduction 

Hybrid systems are widely used in automatic control theory to cope with situations arising when a finite- 
state machine is coupled to mechanisms that can be modeled by differential equations lfT3l . It is the case 
of robots, plant controllers, computer disk drives, automated highway systems, flight control, etc. The 
general behavior of such systems is to pass from one type of smooth dynamics (mode) described by one 
set of differential equations to another smooth dynamics (mode) described by another set of differential 
equations. The command of the modes can be performed by changing one or several discrete variables. 
The mode change can be accompanied or not by jumps (discontinuities) of the trajectories. Depending on 
how the discrete variables are changed, there may be several types of hybrid systems: switched systems 
ll22l . multivalued differential automata ll26l . piecewise smooth systems [91] . Notice that in the last case, 
the mode changes when the trajectory attains some smooth manifolds. In these examples, the changes of 
discrete variables and the evolution of continuous variables are deterministic. The class of hybrid systems 
can be extended by considering stochastic dynamics of both continuous and discrete variables, leading to 
piece-wise deterministic processes, switched diffusions or diffusions with jumps |[T9l l6ll5H23ll4i. Hybrid, 
differential, or stochastic Petri nets provide equivalent descriptions of the dynamics and were also used 
in this context Q. 

The use of hybrid models in systems biology can be justified by the temporal and spatial multi- 
scaleness of biological processes, and by the need to combine qualitative and quantitative approaches to 
study dynamics of cell regulatory networks. Furthermore, hybrid modelling offers a good compromise 
between realistic description of mechanisms of regulation and possibility of testing the model in terms 
of state reachability and temporal logics |[T2l [141 . Threshold dynamics of gene regulatory networks 
El [20l or of excitable signaling systems l28l has been modelled by piecewise-linear and piecewise- 
affine models. These models have relatively simple structure and can, in certain cases, be identified from 
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data |[T8l l8l.Some methods were proposed for computing the set of reachable states of piecewise affine 
models 0. 

Among the applications of hybrid modeling, one of the most important is the cell cycle regulation. 
The machinery of the cell cycle, leading to cell division and proliferation, combines slow growth, spatio- 
temporal re-organisation of the cell, and rapid changes of regulatory proteins concentrations induced by 
post-translational modifications. The advancement through the cell cycle is a well defined sequence of 
stages, separated by checkpoint transitions. This justifies hybrid modelling approaches, such as Tyson's 
hybrid model of the mammalian cell cycle |[24l . This model is based on a Boolean automaton whose 
discrete transitions trigger changes of kinetic parameters in a set of ODEs. The model has been used to 
reproduce flow cytometry data. Instead of building the hybrid model from scratch, another strategy is 
to identify hybrid models from experimental or artificial time series lTT6l fT7l ITTl . The resulting cell cycle 
hybrid models can be used for hypothesis testing, as they are or as parts of larger, integrated models. 

In this paper we develop ideas first introduced in lfl5l . We discuss how a given model of the cell 
cycle, based on ODEs, can be hybridized. The hybridization, based on a tropical geometry heuristics, 
unravels commonalities of cell cycle models. These systems combine quasi-equilibrium states, repre- 
sented by slow invariant manifolds and excitability, represented by rapid transitions to and from these 
manifolds. With respect to fT31 we introduce several new concepts and provide rigorous justification of 
the procedures. Two general hybridization procedures called tropicalizations are introduced in section 2. 
The tropicalized dynamics is guaranteed to be a good approximation for polynomial or rational systems 
with well separated terms and that satisfy a condition called permanency. In subsection 2.2 we introduce 
the tropical equilibration as a method to test permanency. In section 3 we apply these methods to a cell 
cycle biochemical network model. 



2 Tropical geometry and hybridization 

2.1 General settings 

In chemical kinetics, the reagent concentrations satisfy ordinary differential equations: 

dx 

-£=Fi(x),l<i<n. (1) 
Rather generally, the rates are rational functions of the concentrations and read 

F i (x)=P i (x)/Qt(x), (2) 

where P{(x) =EceeA i a /,a a;Q: > Qi( x ) = L/3 e s,A\j3^> ^ e multivariate polynomials. Here a:" = x" l x" 2 ■ ■■x% n , 

x$ = x^'x^ ■ ■■Xn", «j,ce, ^i,B> are nonzero real numbers, and A,-,/?,- are finite subsets of N". 
Special case are represented by 

F i (x)=P+(x)-P i r(x), (3) 

where P^(x), P^(x) are Laurent polynomials with positive coefficients, Pj(x) = HaeA ± a fa x<X ' a ta > 
0, Af are finite subsets of Z". 

In multiscale biochemical systems, the various monomials of the Laurent polynomials have different 
orders, and at a given time, there is only one or a few dominating terms. Therefore, it could make sense 
to replace Laurent polynomials with positive real coefficients Y,aeA a a& a , by max-plus polynomials 
max a€A {log(a a )+ < log(x),a >}. 
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This heuristic can be used to associate a piecewise-smooth hybrid model to the system of rational 
ODEs (Q]), in two different ways. 

The first method was proposed in |[T5l and can be applied to any rational ODE system defined by 
©,©: 

Definition 2.1 We call complete tropicalization of the smooth ODE system (Q]), (O the following piecewise- 
smooth system: 

dx- 

-j = DomPi(x) /DomQi(x) , (4) 

where Dom{a i: aX a } aeA , = sign{ai, amax )exp[max a eA,{ lo s{\ a i..a\)+ <u,a >}}. u = (logx u . . . ,logx n ), 
and ai^a max , CC ma x £ M denote the coefficient of the monomial for which the maximum is attained. In simple 
words, Dom renders the monomial of largest absolute value, with its sign. 

The second method,proposed in ETI . applies to the systems CD,©. 

Definition 2.2 We call two terms tropicalization of the smooth ODE system (Q~|), (J3J) the following piecewise- 
smooth system: 

dx- 

— = DomP? (x) - DomPr (x), (5) 
dt iv/' 

The two-terms tropicalization was used in II2T1 to analyse the dependence of steady states on the model 
parameters. The complete tropicalization was used for the study of the model dynamics and for the 
model reduction |fT5i 

For both tropicalization methods, for each occurrence of the Dom operator, one can introduce a 
tropical manifold, defined as the subset of R" where the maximum in Dom is attained by at least two 
terms. For instance, for n = 2, such tropical manifold is made of points, segments connecting these 
points, and half-lines. The tropical manifolds in such an arrangement decompose the space into sectors, 
inside which one monomial dominates all the others in the definition of the reagent rates. The study of 
this arrangement give hints on the possible steady states and attractors, as well as on their bifurcations. 

2.2 Justification of the tropicalization and some estimates 

In the general case, the tropicalization heuristic is difficult to justify by rigorous estimates, however, this 
is possible in some cases. We state here some results in this direction. To simplify, let us consider the 
class of polynomial systems, corresponding to mass action law chemical kinetics: 

dx- M 

= Pi(x,e) = £ Mjj(x, e), Mij = ^(e)x^ (6) 

where a, ; are multi-indices, and £ is a small parameter. So, the right hand side of © is a sum of 
monomials. We suppose that coefficients Pjj have different orders in e: 

P i j(e) = s b -Jp u , (7) 

where h i} / b fj > for (ij) / (/',/). 

We also suppose that the cone R> = {x : x,- > 0} is invariant under dynamics © and initial data are 
positive: 

Xi (0) >S>0. 

The terms (|7]) can have different signs, the ones with Pjj > are production terms, and those with < 
are degradation terms. 
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From the biochemical point of view, the choice (0 is justified by the fact that biochemical processes 
have many, well separated timescales. Furthermore, we are interested in biochemical circuits that can 
function in a stable way More precisely, we use the permanency concept, borrowed from ecology (the 
Lotka -Volterra model, see for instance 11251 ). 

Definition 2.3 The system © is permanent, if there are two constants C > and C + > 0, and a function 
To, such that 

C < Xj(t) < C + , for all t > Tq(x(0)) and for every i. (8) 
We assume that C± and Tq are uniform in (do not depend on) £ as £ — > 0. 

For permanent systems, we can obtain some results justifying the two procedures of tropicalization. 

Proposition 2.4 Assume that system ([6]) is permanent. Let x, x be the solutions to the Cauchy problem 
for ((6]) and (0) (or (O ), respectively, with the same initial data: 

jc(0) = x(0). 

Then the difference y(t) = x(t) —x(t) satisfies the estimate 

\y(t)\<C ie y e x V (bt), 7 >0, (9) 

where the positive constants C\ , b are uniform in £. If the original system ((6]) is structurally stable in 
the domain Qc-,c + = { x '■ C < \x\ < C + }, then the corresponding tropical systems (0) and (0 are also 
permanent and there is an orbital topological equivalence x = h e (x) between the trajectories x(t) and 
x(t) of the corresponding Cauchy problems. The homeomorphism h e is close to the identity as £ — > 0. 

The proof of the estimate © follows immediately by the Gronwall lemma. The second assertion follows 
directly from the definition of structural stability. 

Permanency property is not easy to check. In the case of systems © we can make a renormalization 

Xi = £ a 'Xi (10) 

and suppose that ([8]) holds for Xj with Cf uniform in e. 

We seek for renormalization exponents a, such that only a few terms dominate all the others, for 
each i-th equation © as e — > 0. Let us denote the number of terms with minimum degree in e for i-th 
equation as ra,-. Naturally, 1 < m,- < M,-. After renormalization, we remove all small terms that have 
smaller orders in e as £ — > 0. We can call this procedure tropical removing. The system obtained can be 
named tropically truncated system. 

Let us denote a'/ the I th coefficient of the multi-index dij. If all m; = 1 then we have the following 
truncated system 

^ = e**Kt), ^= A7(0 x« y(i) , (11) 

where 

Mi = Yij(f) L ai <<H (12) 

and 

n 

Mi >YijY, a l Jai forall j^j(i). (13) 
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If all nij = 2, in order to find possible renormalization exponents a,-, it is necessary to resolve a family of 
lineal - programming problem. Each problem is defined by a set of pairs (j(i),k(i)) such that j(i) ^ k(i). 
We define jU, by 

n ■ -i -\ n 
IM = Yij(i) + £ a\ ai = 7ik(f) + £ a i a i ( 14 ) 
i=\ 1=1 

and obtain the system of the following inequalities 

n 

IM > Yij + £ a'/m for all j ^ ,k(i) . (15) 
i=i 

The following straightforward lemma gives a necessary condition of permanency of the system ©. 

Lemma 2.5 Assume a tropically truncated system is permanent. Then, for each i £ {1, . . . ,n}, the i-th 
equation of this system contains at least two terms. The terms should have different signs for coefficients 
Pip i.e., one term should he a production one, while another term should he a degradation term. 

We call "tropical equilibration", the condition in Lemma [231 This condition means that permanency is 
acquired only if at least two terms of different signs have the maximal order, for each equation of the 
system ©. This idea is not new, and can be traced back to Newton. 

The tropical equilibration condition can be used to determine the renormalization exponents, by the 
following algorithm. 

Step 1. For each i let us choose a pair (j(i),k(i)) such that {1, . . . ,ffl,-} and j < k. The sign of 
the corresponding terms should be different. 

Step 2. We resolve the linear system of algebraic equations 

n . n 

Yij(i) ~ Yik(i) = - £ ot'/ {,) a, + £ of {l) ai, (16) 
i=i 1=1 

for ai, together with the inequalities (fT3T) . 

Notice that although that Step 2 has polynomial complexity, the tropical equilibration problem has 
an exponential number of choices at Step 1. 

Assume that, as a result of this procedure, we obtain the system 

d Ji = s»(F+(x)-F i -(x)), Ff=pt^. (17) 
One can expect that, in a "generic" casd3, all /i; are mutually different, namely 

= Ml <jU2 < ••• < Hn-l < Hn- (18) 

We can now state a sufficient condition for permanency. Let us consider the first equation (fTTT ) with i = 1 
and let us denote y = x\,z = (x2,...,x n ) tr . In this notation, the first equation becomes 

^=f(y)=b l (z)y pi -b 2 (z)y fh , b h b 2 >0, ftsR. (19) 

Since \i 2 > 0, one has that z{t) is a slow function of time and thus we can suppose that bj are constants 
(this step can be rendered rigorous by using the concept of invariant manifold and methods from ifTTTl ). 

1 supposing that multi-indices cty are chosen uniformly, by generic we understand almost always except for cases of vanish- 
ing probability, see also 1101 
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The permanency property can be then checked in an elementary way. All rest points of ( fT9l ) are roots of 
/. If / > 0, y(t) is an increasing time function and if / < 0, y(t) is a decreasing time function. A single 
root y\ of / within (0, +oo) is given by 

b\ Pi - p 2 

These properties entail the following result: 

Lemma 2.6 Equation \19\has the permanence property if and only if 

0<ft</3 2 , i 

or 

ft < ft < 0, ii 

or 

ft < 0, ft > 0. iii 

For fixed z, in these cases we have 

y(t,z)^yo, ast^oo. 



The generic situation described by the conditions (PT8T ) lead to trivial "chain-like" relaxation towards a 
point attractor, provided that we have permanency at each step. This result is the nonlinear analogue of 
the similar result that monomolecular networks with total separation relax as chains and can only have 
stable point attractors ifTOll . 

The following theorem describes a less trivial situation, when limit cycles are possible. 

Theorem 2.7 Assume = pi\ < \l 2 < ... < < ]X n holds. If the procedure, described above, leads to 
the permanency property at each step, where i= 1 , 2, — 2, and the last two equations have a globally 
attracting hyperbolic rest point or globally attracting hyperbolic limit cycle, then the tropically truncated 
system is permanent and has an attractor of the same type. Moreover, for sufficiently small e the initial 
system also is permanent for initial data from some appropriate domain W em ,a an d nas an analogous 
attracting hyperbolic rest point (limit cycle) close to the attractor of the truncated system. If the rest 
point (cycle) is not globally attracting, then we can say nothing on permanency but, for sufficiently small 
E, the initial system still has an analogous attracting hyperbolic rest point (limit cycle) close to the 
attractor of truncated system and the same topological structure. 

Finally, let us note that tropical equilibrations with permanency imply the existence of invariant mani- 
folds. This allows to reduce the number of variables of the model while preserving good accuracy in the 
description of the dynamics. The following Lemma is useful in this aspect. 

Lemma 2.8 Consider the system 

^=f(y)=bi(z)y pi -b 2 (z)y^, b u b 2 >0, /3, e R. (21) 



(22) 
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where z € R m , A > is a parameter and the function F enjoys the following properties. This function 
lies in an Holder class 

F G C 1+r , r > 0, 

and the corresponding norms are uniformly bounded in Q = (0, +°°) x W, for some open domain W C 
R m 

\ F \c l + r (Q) < C 2- 

Assume one of conditions i, ii, iii of Lemma \2.6\ holds. We also suppose that hi are smooth functions of z 
for all z. such that \z\ > 5q > 0. Assume that z. € W implies \z,\ > 8q. 

Then, for sufficiently small X < Xo(C 2 ,bi,b2,fii,fh.,8) equations rt2il ), A23\) have a locally invariant 
and locally attracting manifold 

y = Y(z,X), Y€C 1+r (W), (23) 

and Y has the asymptotics 

Y(z,X)= yi (z)+Y, Y€C 1+r (W), (24) 

where 

\f(z,X)\ c i+r {w) <C S X S , s>0. (25) 



Proof. The proof is standard, follows from Theorems in ifTTTl . Ch. 9. 

3 A paradigmatic cell cycle model and its tropicalization 

3.1 Description of the model 

We study here the cell cycle model proposed by Tyson ll27l . This model mimics the interplay between 
cyclin and cyclin dependent kinase cdc2 (forming the maturation promoting factor MPF complex) during 
the progression of the cell cycle. The model demonstrates that this biochemical system can function as an 
oscillator, or converge to a steady state with large MPF concentration, or behave as an excitable switch. 
The three regimes can be associated to early embryos rapid division, metaphase arrest of unfertilized 
eggs, and growth controlled division of somatic cells, respectively. This model takes into account auto- 
catalytic activity of MPF (positive feedback). It can be described as a nonlinear cycle of biochemical 
reactions and corresponds to the following set of differential equations: 

y'\ = e~ 3 k 9 y 2 - e~%yi + k 6 y 3 , y' 2 = e^ 6 k s yi - s~ 3 k 9 y 2 - S~ 2 k 3 y 2 y5, 

y' 3 = S 2 k' 4 y4. + e^ 2 k 4 y 4 yl - k 6 y 3 , y' 4 = -S 2 k' 4 y 4 - E~ 2 k 4 y 4 y\ + £ _2 ^ 3 3 ; 23 ; 5 , 

y' 5 =e 2 k l -e- 2 k 3 y 2 y 5 , (26) 

Here k t > are rate constants, y j: i E [1,5] are concentrations of cdc2, p-cdc2 (phosphorylated kinase), 
cyclin-p:cdc2 complex (active MPF), cyclin-p:cdc2-p complex (inactive MPF), and cyclin, respectively. 
With respect to the original model we have introduced a small parameter e > to cope with the order of 
the rate constants (e = 0.1 in the original model). 
The system (l26l ) has the conservation law 

yi(t)+y2(t)+y 3 (t)+y 4 (t) = l, (27) 

where the value 1 (total initial concentration of kinase cdc2) was chosen by convenience. 
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3.2 Tropical equilibrations and model reduction 

Let us apply the tropical equilibration principle. To this aim, we renormalize the variables, 

yi = e a %. (28) 
Let us substitute these relations into the system of equations. As a result, we obtain 

y\ = £-i+«2-« lk9 y 2 _ e -6 km + k 6 e a ^y^ y' 2 = e^^hm - e~Xn - £- 2+a5 k 3 y 2 y 5 , 

y 5 =e 2 -^k l -e- 2+ ^k 3 y 2 y 5 . (29) 

The system (1291 has the conservation law 

e fll y 1 (f) + e*y 2 (/) + 6^(0 + 6^4(0 = 1. (30) 

In order to compute the exponents a, we use tropical equilibrations together with the conservation law 
(|3Ql >. There are 2 4 valiants of tropical equilibrations. To our surprise, there is only one solution for the 
exponents values. We can show that all possible equilibrations of the variables y 3 , y$ and y$ uniquely set 
the values of two exponents, a 3 = 2,a^ = 0. 

Let us consider the variants with respect to the equilibrations of the variables y\ and y 2 . Denoting by 
Tj the i th term in the equation, we have the following situations: 

1) Ineq. for ft: T\ = T2, T3 <= Tl, In eq. for y 2 : Tl = T2, T3<=T\. 

2) In eq. for ft: T\=T2, T3 <= Tl, In eq. for y 2 : T1 = T3, T2<=T3. 

3) In eq. for ft: T2 = T3, Tl <= T2, In eq. for ft: T1 = T2, T3<=T1. 

4) Ineq. for ft: T2 = T3, Tl <=T2, In eq. for ft: Tl = 73, T2<=T1. 

In the variant 1 (Case I) the tropical equilibrations do not fix the values of the exponent as and we 

get 

Cl\=l — <Zg, £?2 = 4 — fl5, (33 = 2, C?4 = 0, £?5 > — 1. (31) 

In the variants 2,3,4 (Case II) the exponents are uniquely determined from equilibrations and we obtain 

a\ = 8, a 2 = 5, a 3 = 2, a\ = 0, a$ = — 1. (32) 

However, (l32l and the conservation law (l30l are incompatible, therefore, the case II can be rejected. 
In the case I the conservation law takes the form £ 7 ~" 5 ft (t) + £ 4 ~ fl5 ft(f) +ft(f) = 1 + o(l), as £ — > 0. 
Assuming that ft (0) = 0(1) and ft ^ 1 (it is reasonable, since it is a "generic case"), we obtain = 4. 
Thus, the only possible situation is variant 1 (Case I) and the corresponding set of exponents is: 

a\ = 3, a2 = 0, £?3 = 2, c?4 = 0, <35 = 4. (33) 

Let us note that the terms Tl and T2 in the equations for the variables y\,y 2 correspond to direct and 
reverse rates of a phosphorylation/dephosphorylation cycle transforming y\ into y 2 and back. Thus, 
biochemically, (Case I) corresponds to the quasi-equilibrium of this cycle. Furthermore, the equilibration 
of all the variables leads to the exponents (f3TT >. In this case, 2 + 04 — a 3 = —2 + a 3 + 04 = 0, 2 = 
—2 + 2fl3 = —2 + a 2 + a=, — a\, meaning that the tropical equilibrations of the variables y 3 , J4 are triple 
(in each equation, all three terms have the same order). 
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We finally obtain the following renormalized system 

y\ = e- 6 (k 9 y 2 -hy 1 )+k 6 £- 1 y 3 , y' 2 = e~ 3 {hyi -hn) -e 2 %Ms. 
y 3 = k' 4 y 4 + k 4 y 4 yl~k 6 y 3 , y 4 = E 2 (-k' 4 y 4 -&4M3 + k 3 y 2 ys), 
y' 5 =e- 2 (k l -k 3 y 2 y 5 ). (34) 

The structure of the system 04l emphasizes the multiple time scales of the model. The fastest variables 
are in order y\, then y 2 and y 5 . The variables y 3 and y 4 are slow. 
Assume that 

y 2 > 5 > 0. (35) 

This important assumption ensures the existence of an invariant manifold and will be justified, a posteri- 
ori. 

Then, from the last equation 04l we obtain the relation 

y 2 y5=k l /k 3 + 0(e 2 ), 

which represents the equation of an invariant manifold. 

In turn, this relation leads to the following equations for y 3 ,y 4 

y' 3 =k' 4 y 4 + k 4 y 4 y 2 -k 6 y 3 , y' 4 = £ 2 {-k' 4 y 4 - k 4 y 4 y\ + k\). (36) 

A second invariant manifold equation is defined by the equation 

yi=k^(k 9 y 2 +k 6 e 5 y 3 ). (37) 

Remind that 3^3 is a slow variable. Then, for y 2 we have 

y' 2 = e 2 {hy 3 -h). (38) 

System d36l ) represents a two-dimensional reduced model of the initial five-dimensional system. This 
result shows that tropical equilibrations can be used for model reduction. 
The solutions of (l36l) either tend to the stable equilibrium 

k' 4 + k 4 {ki/k 6 y 

or, if this equilibrium is unstable, to a limit cycle. 

Based on the general Theorem [277] we can assert the following: 

Theorem 3.1 Assume ( T33T ) holds with 8 > 0. If the shorted system 06l) has a stable hyperbolic limit cy- 
cle, then, under above conditions, for sufficiently small e the five component system d26D also has a stable 
limit cycle. If the shorted system (136\) has a stable hyperbolic equilibrium, then, under above conditions, 
for sufficiently small E the five component system \26\ also has a stable hyperbolic equilibrium. 



We have studied the system 061) analytically and numerically. The numerical simulations confirm the 
criteria of cycle existence both for small e and for e = 0(1). For small epsilon the cycle has a singular 
structure. The amplitude of y 3 and the cycle period increase in e, approximative^, as e~ 2 (the assertion 
about period is natural since the rate of y 4 is 0{e 2 )). 
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Hyperbolicity can be straightforwardly checked for the rest point (|39l , by computing the eigenvalues 
of the linearized system. Denote by Y = (y$, y 4 ). When the rest point Y° = (^3,^4) is hyperbolic and 
stable we have the following estimate 

\Y(t)-Y°\ <Ciexp(-ci£ 2 ?) (40) 

with some C\ , c\ > holds. Integrating (l38l) for yi over interval [0, t] gives 

|y 2 (T)-y 2 (0)| <e 2 C lC - 1 =o(l) (41) 

uniformly in % > as e — > 0. This yields that y2(t) > 8 if ^2(0) > 28 and therefore, y 2 (0 does not go to 
zero for large t, justifying the estimate (I35T ) needed for the existence of an invariant manifold. The case 
of a limit cycle is discussed in the next subsection. 

3.3 Singular limit cycle and hybrid dynamics 

Up to this point, the tropical ideas were used for reducing the dynamics of the model. In this section we 
show that the tropicalization heuristic is well adapted for decomposing the limit cycle into slow and fast 
modes, providing a hybrid description of the dynamics. 

Let us note that in a hybrid, excitable system, it is possible that not all variables are equilibrated. 
Also, the system can have more than two different equilibrations and associated invariant manifolds, and 
jump from one invariant manifold to another during the dynamics. Let us consider that the variables 
yi,y2,y5 are equilibrated as above, but now, only one among the variables 373 or y 4 are equilibrated. We 
have four situations: 

l)Ineq. for ^3: Tl = T3, T2 <= Tl, 2) In eq. fory 3 : T2 = T3, Tl <=T2, 

3)Ineq. fory 4 : Tl = T3, T2 <= Tl, 4) In eq. fory 4 : T2 = T3, Tl <=T2. 

Combined with the conservation law condition (l30l) . variants 1 and 2 lead to the same triple tropical 
equilibration as before (Case I) and exponents 03]). We denote the corresponding invariant manifold 
^#1. The renormalized equations are the same as (l36l ). 

Variant 3 can be rejected by the general permanency criterion given by Lemma [2T61 Variant 4 (Case 
III) satisfies the permanency criterion and leads to 

ai=3, t?2 = 0, £23 = 0, ^4 = 4, =4. (42) 

This corresponds to a double equilibration (two equal terms) of the variable y\, the variable y^ being not 
equilibrated. We denote the corresponding invariant manifold ^#2- The renormalized equations read 

y 3 =£ 6 k' 4 y 4 + £ 2 k4y 4 y]-k6h, ^ = -e^ + r^-^M^+^i)- (43) 

We can provide a hybrid description of the cell cycle, by decomposing the periodic orbit into three 
modes (FigQ]). The first mode is the slowest and has the longest duration. It consists in the dynamics on 
the slow invariant manifold ^#1 at low values of 3^3, and can be described by the algebraic-differential 
system y' 3 = £ 2 &i — ^6^3 , ^4^4 + fco^J 2 — = (part between 0\ and O of the orbit in FigQ]:)). The exit 
from the invariant manifold -#1 occurs at a critical point (point O). The next slowest mode conesponds 
to the decrease of yi (part between O2 and 0\ of the orbit in FigQ};)) and can be described by two terms 
truncated system y' 4 = £~ 2 {k\ — k^y^y^), 3^ = — Finally, there is a fast mode, corresponding to the 
fast increase of y'3 (part between O and O2 of the orbit in FigQ]:)) and described by the truncated system 
y' A = —£ 2 k 4 yny 2 , y' 3 = e~ 2 k 4 y 4 y 2 . One can notice (FigQ]:)) that this hybrid approximation is very accurate 
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for small e. At a distance from the tropical manifold, the hybrid orbit coincides with the one generated 
by the two terms or by the complete tropicalization. However, close to the tropical manifold, the two 
term and the complete tropicalization are less accurate than the hybrid approximation described above. 

Below we state rigorous estimates describing the slow movement on ^#1 and the fast jump towards 
^#2- The two terms description of the dynamics on is a direct consequence of the Proposition 12.41 
and Lemmas 12.612.81 

To simplify notation, we rewrite the system (l36l) for y-$-,y\ as 

x =y + yx 2 -k x, (44) 
/ = e 2 (-y-yx 2 + k l ). (45) 

We can obtain such a presentation by a linear variable change. 
Let us define the functions 



K ' 2y 



k 2 + Jk 2 -4y 2 

and the points 

y = k /2, x =X(y ) = l. 

Lemma 3.2 Solutions of fi44\. ri?5l ) with initial data x(0),y(0) such that 

0<$)<x(0) <X+(y )-So, 0<y(0)<y -8 , (46) 
where 8q is a small positive number independent on £, satisfy 

\x{t)-X{y{t))\<Cy{£ + exv{-c y 8 t)), t < T (x(0),y(0),e) (47) 
this estimate holds while 

y(t)<yo-S2, o\ >0. (48) 
Let us find some estimates of solutions at y = yo,x = xq. Our goal is to prove that at this point x(t) starts 



to increase sharply. After this, the terms iyx 2 play the main role in the equations (|441.(|43T). and the other 
terms can be removed while the x-component is big. 
Let us introduce new variables u,v by 



x-x = v, u = y-y Q . 

For u,v one obtains 

v' = ^v 2 + M (l + (l+v) 2 ), (49) 



u = £ 2 (-(| +«)(1 + (1 +v) 2 ) + h)) = e 2 g(u,v). (50) 
Let us consider for this system the Cauchy problem with initial data 

v(0)=v , «(0)=«o. (51) 
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Lemma 3.3 Consider the Cauchy problem < \49\) , ( |50| ) and ([37]) under assumptions that 

u = k>0, \v \ < S 5 , 

and 

ki>k , (52) 

where 80,85 are small enough (but independent on e). Let A be a large positive number independent of 
e. Then within some interval t S (^0(85, ko,ki,A),Ti(S5,ko,ki,A)) one has 

u(t)>0, v(t)>A, v(z )=A. (53) 

This result can be reinforced. Actually, x(t) attains values of the order 0(e~ 2 ). 
Lemma 3.4 Assume 

x(h)=A»l, y(h)>ko/2. (54) 

Then 

x{t)>{A- l - l -at)-\ (55) 
y(t) > a, (56) 

for <r,A such that 

oA > 2k , a > fc exp(-2a _1 ), (57) 



and t such that 



4 Conclusion 



(A 1 - |cjO > e 2 . (58) 



We showed that tropical ideas can be usefully employed to reduce and hybridize polynomial or rational 
dynamical systems occurring in modelling the molecular machinery of the cell cycle. The main idea 
consists in keeping only the dominant monomial terms in the right hand side of the ordinary differential 
equations. Depending on the position in phase space, one should keep one, two, or more such terms. 
The places where two or more monomial terms are equal form the so-called tropical manifolds. The one 
term approximation is valid far from the tropical manifolds, whereas close to tropical manifolds several 
dominating terms of opposite signs can equilibrate each other. These "tropical equilibrations" of the 
dominating terms slow down the dynamics and produce attractive invariant manifolds. 

The possible applications of this method are multiple. Generally, the method can be used to obtain 
simplified models. In the example studied here, we have started with a five variables model, that has been 
reduced to two variables and hybridized. The modes of the hybrid model have the simple structure of 
monomial differential or differential-algebraic equations. Two general methods that we called complete 
and two terms tropicalizations provide description of the modes and of the mode changes. However, 
these general procedures may lead to inaccurate approximations when the full model does not satisfy 
permanency globally. In such cases, more thorough analysis is needed. We have shown that the model of 
embryonic cell cycle has essentially three modes with different timescales, namely slow accumulation of 
cyclin, rapid activation of MPF and intermediately rapid degradation of cyclin and inactivation of MPF. 
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Figure 1: Limit cycle behavior of the paradigmatic cell cycle model from ll27l . (a) The three main 
processes during the embryonic cell cycle are, in order of the timescales, fast increase of yj (active MPF, 
triggering mitosis), slower decrease of V3 and very slow increase of y\ (inactive MPF). (b) Two invariant 
manifolds corresponding to the two slow processes are close to the tropical manifolds (blue lines) and 
result from equilibration of the variables (equilibration of y^ corresponds to ^2 and equilibration of 
corresponds to (c) A three modes hybrid approximation of the cell cycle (in red) compared to the 
original limit cycle (black crosses). 
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The fastest mode is described by monomial ODEs, whereas the less fast modes correspond to tropical 
equilibrations and are described by differential-algebraic equations. 

Several improvements and developments are needed in order to apply these methods at a larger 
scale. The computation of tropical equilibrations suffers from combinatorial explosion. However, for the 
biochemical network used as working example, the number of solutions seems to be very small compared 
to the large combinatorics of monomial terms. There is hope, that once formulated in constraint logic 
programming, the problem of equilibrations could be efficiently computed in practice as a constraint 
satisfaction problem. Also, effective methods are needed to compute the transitions between modes. 
The main difficulty here is related to walls (segments of the tropical manifolds) crossing. Near walls, 
two or more terms are dominant. When these terms are equilibrated, orbits remain close to the walls 
and are contained in invariant manifolds. The complete or two terms tropicalizations provide general 
heuristic for mode transitions. These approximations mail fail close to walls. For instance, as we showed 
in a previous paper [fT5l , the complete tropicalization predicts sliding modes that evolve on the wall and 
stay thus close to orbits of the full system. However, these sliding modes can be too long, leaving the 
wall when the orbits of the full system are already far away. In order to get an accurate description of 
the behavior near such walls we had to compute invariant manifolds. Although this is generally much 
simpler than integrating the full set of equations, it could become difficult for tropical equilibrations 
involving more than two terms. Future work will be dedicated to developing general methods for this 
problem. 
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5 Appendix: proofs 



Proof of Lemma 3.2 Let us consider the equation 

x =y(0)+y(0)x 2 -k x = f (x). 



(59) 
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One observes that 

/(*) < S 4 < 0, x G (X(y(0)) + 8o,X + (y ) - do), 
f{x)>8 4 , xe (8o,X{y )-8o), 

and 

f'(x)\ x=x{yo) <-8 5 <0. 

Therefore, if < x(0) < ^ 0+ v^o^ 4 y W_ _ ^ t j ien ^ so i u ti n x(t) attains a small 5- neighborhood of 
X(y(0)) within a bounded time interval T\ (8, So, 84): 

|*(/,*,,yd)-XCy(0))| < 8, t = Z\ (5, 

Within a small 8 - neighborhood of X(y(0) we set u = x — X(y(0)) and then we can rewrite (l59l ) as 
follows : 

u' = -Ku + h(u), |/i(m)|<Cim 2 , k(7i) = 5. (60) 
where K" > is independent of 5. Then, if 8 is small enough, we have that u{t) < 8 for all t > T\ and 

KOI < C(7b)exp(-Kf/2), t > To. 

Let us compare now the solution x{t) of d59l ) and the coiTesponding solution x(t,xo,yo) of (l44l) with the 
same initial data. For x(t) — x(0) = w one has, since y(t) —y(0) < C{To)e 2 on any bounded interval 

te[o,T ], 

w t =a(t)w + e 2 g(x,t,w),0 <t <T Q . 
with a smooth function g, a(t) is bounded function. Now the Gronwall inequality implies 

K*)| <C 2 (T Q )e 2 , t€[0,T ]. 

Therefore, one has 

\x(t,XQ,y )-X(y(t))\ <8, t = T 2 (8, 80,84). 

Since the function X (y) defines a smooth, locally attracting (for y < yo) invariant manifold, this proves 
our assertion. One can prove this in another, elementary way. Let us define u = x — X(y). Then 

u' = -K(t)u+h(u) + 0(e), \h(u)\<Ciu 2 , u(Ti) = 8, (61) 

where rc(t) > Kq while y < y — So- Again one has u(t) < 28 for t > T\ (while y < yo — So). Thus, 

u' < - yW + 0(£) 

that entails the need estimate |47] 

Proof of Lemma 3.3 Let us consider the Cauchy problem 

w/ = ^V, w (0) = v = v(0) (62) 
It is clear, by the comparison principle, that 

v(t) > w(t) 
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while u(t) > 0. Consequently, the assertion vv(Tb) = A proves the lemma. Let us prove first that if 8(, is 
small enough (but independent on e), for some t = X\ we have 

v(Ti) = fc. (63) 

Without loss of generality, we assume that 

v(t)<8 6 , 0<*<Ti, v(ti) = 5 6 . (64) 
Let 8& > be small enough such that 

r "(l + (l + fc) 2 )+*i>0. (65) 



2 

Such a choice of S(, is possible due to d52l ). 

Assume u(t) > within some time interval [0, T2], and ufa) — f° r some T2 < T\. Since w(0) = 0, 
we have T2 > 0. We can suppose without loss of generality that T2 is the first moment, where ufa) = 0. 
Then 

u'(x 2 )<0. 

But then we obtain a contradiction with (|65T ) at ? = T2, since the right hand side of this equation is positive 
at this time moment. 

Therefore, we have shown that u(t) > for all t from [0, Ti] if v(t) < 8e for such t. Then v(t) > w(t) 
on this time interval. The function w can be found, and an easy computation gives 

w (t) = ( W -\0)-2k [ ty l ) = {u -2k l t). (66) 

Assume that v(t) < 8(, for all t. Then (l66l ) holds for all t, but w(Z\) > 8(, for some Z\. We have obtained 
a contradiction, thus (l63l) is proved. 

Let us prove that v(t) = A. Let us consider an interval [t\ , T] such that 

|v(0|<A, ?G[?i,r], v(h) = 8 6 . (67) 

Then ( f65l > implies 

« / <£ 2 ((y + |«|)(l + (l+A) 2 )+^ 1 )), (68) 
that gives, by the Gronwall lemma, 

u(t) <C 2 u(0)exp(Ci(A)e 2 t). (69) 

Within interval (t\ , T] one has then 

u{t)<C 4 K, te[0,T]. (70) 

and, therefore, 

v > yV 2 -C 4 fC, t€[0,T]. (71) 

Suppose that 

C4K < —06- 



Then dH) entails 

v'> 

4 



v'>^v, ?G[0,r], (72) 
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This gives 

v(f)>fcexp(^(f-fi)). 

This leads to a contradiction for T large enough (however, let us remark that T is uniform in e). 
Proof of Lemma 3.4 

Remark: to satisfy (l57l ) it suffices to set a = 2koA with a large A. 

Suppose that either the estimate (|55T ) (case A) or the second estimate d56l ) (case B) is violated at some 
T but the both inequalities hold for all t\ < t < T. We can set t\ = 0. Let us consider the case A. Since 
(1561 ) hold, one has 

x>ox 2 , x(0)=A, t€[0,T]. 

This implies 

x'/x 2 = (-l/x)' > o, x(0) =A. 

Thus, 

x> (A^-at)' 1 . 

For t = T this result gives a contradiction with (T55T ). 
Let us consider the case B. Since (|55"I ) hold, one has 

y\t)<-£ z (y{A- l -ot)- x -k x ), y(0)=k Q /2, te[0,T]. 

This implies 

y(t)< ^exp(-£ 2 f (l + (A- 1 -os)~ 2 )ds) + e l ki Texp(-£ 2 f \l + (A~ l -as)" 2 )ds)dx. 

2 Jo JO Jt 

for t € [0,T]. Notice 

/ T f (l + (A" 1 - os)- 2 )ds = t -x + o-\{A- Y - at)- 1 - (A" 1 - at)- 1 ). 
Jo 

Thus, for sufficiently small e, taking into account (T58T ) one has y(T) < -^-(-2/a). Under condition (E7J 
this result gives a contradiction with d56l ). 
Proof of Theorem 2.7 

i Suppose that the tropically truncated system (TTS) has a globally attracting compact invariant set 
si ' . Let be JT be an open neighborhood of this set. We can choose this neighborhood as a box that contains 
si '. Then, for all initial data x(0), the corresponding trajectory x(t),x(0) lies in JT for all t > Tq(xq, JT). 
Therefore, our TTS is permanent. Here we do not use that the cycle (rest point) is hyperbolic. 

ii Permanency of the initial system follows from hyperbolicity of s/. Hyperbolic sets are persistent 
(structurally stable) (Ruelle 1989). Since this set is globally attracting, all TTS is structurally stable (as a 
dynamical system). This implies that the initial system has a hyperbolic attractor close to s/, since initial 
system is a small perturbation of TTC in JT. 

iii If the set s/ is only locally attracting, the last assertion of Theorem follows from persistency of 
hyperbolic sets. 



